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ABSTRACT 

We present the first Lyman Alpha Emitter (LAE) study that combines: (i) cosmolog- 
\ ical SPH simulations run using GADGET-2, (ii) radiative transfer simulations (CRASH), 

and (iii) a previously developed LAE model. This complete LAE model accounts for 
the intrinsic LAE Lya / continuum luminosity, dust enrichment and Lya transmission 
, <~| through the intergalactic medium (IGM), to quantify the effects of reionization, dust 

Qh and velocity fields on the Lya and UV Luminosity Functions (LF). We find that a 

model neglecting dust sorely fails to reproduce either the slope or the magnitude of 
£^ . the observed Lya and UV LFs. Clumped dust is required to simultaneously fit the 

£^ ' observed UV and Lya LFs, such that the intrinsic Lya-to-continuum luminosity is 

enhanced by a factor f a /f c ~ 1-5 (3.7) excluding (including) peculiar velocities. The 
higher value including velocity fields arises since LAEs reside in large potential wells 
and inflows decrease their Lya transmission. For the first time, a degeneracy is found 
between the the ionization state of the IGM and the clumping of dust inside high- 
redshift galaxies. The Lya LF z ~ 5.7 can be well reproduced (to within a 5a error) 
I by a wide range of IGM average neutral hydrogen fraction, 3.4 x 1CT 3 < (xhi) < 0.16, 

OO ■ provided that the increase in the Lya transmission through a more ionized IGM is 

■ compensated by a decrease in the Lya escape fraction from the galaxy due to dust 
I absorption. The physical properties of LAEs are presented, along with a discussion of 

■ the assumptions adopted. 
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1 INTRODUCTION models as fresh data sets are acquired. In this sense, it 



The epoch of reionization marks the second major change 
in the ionization state of the Universe. Reionization begins 



has been suggested (Malhotra & Rhoads 2004, 2005; San- 
tos 2004; Mesinger, Haiman & Cen 2004; Haiman & Cen 

2005; Dijkstra, Lidz & Wyithe 2007; Mesinger & Furlanetto 
when the first sources of neutral hydrogen (H i ) ionizing 200g . Dayal) perrara & GaUerani 2008; D 1 et aL 2009a) 

photons form withm dark matter potential wefls and start „„.„, x , , r 1 ■ 1 1 1 -a 1 hit 

; iiii 2009 b) that a cfass of high redshift gafaxies cafled Lyman 

bunding an ionization region around themselves, the so- . , , „ , /T , „ s , , , . , 

° ° Alpha Emitters (LAfis) coufd be an important addition to 

caffed Stromgren sphere. However, the reionization history , , , . . . . . , . 

, , compfementary data sets to constrain the reionization his- 

and the redshift at which it ends stiff remain the subject 

of much discussion. This is because the reionization process 

depends on a number of parameters inciuding the initiai LAEs > gafaxies identified by means of their very strong 

mass function (IMF) of the first sources, their star forma- L ya ffne (1216 A) emfssion, have been becoming increasingfy 

tion rates (SFR), their steffar metafficity and age, the escape popufar as probes of reionization for three primary reasons, 

fraction of H 1 ionizing photons produced by each source and Firstiy, speeffic signatures fike the strength, width and the 

the ciumping of the intergaiactfc medium (IGM), to name a continuum break bfuewards of the Lya fine make the detec- 

f ew . tion of LAEs unambiguous to a iarge degree. Secondfy, sfnee 

Given the iarge number of free parameters that in- L J a photons have a large absorption cross-section against 

evitabiy enter into the construction of theoreticai reioniza- H 1 » th eir attenuation can be used to put constraints on the 

tion modeis, it is imperative to compare and update the ionization state of the IGM. Thirdly, there are hundreds of 

confirmed LAEs at z ~ 4.5 (Finkeistein et ai. 2007), z ~ 5.7 
(Malhotra et af. 2005; Shfmasaku et af. 2006) and z ~ 6.6 

* E-mail: dayal@sissa.it (PD) (Taniguchi et ai. 2005; Kashikawa et al. 2006), which ex- 
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actly probe the redshift range around which reionization is 
supposed to have ended. 

The data accumulated on LAEs shows some very sur- 
prising features; while the apparent Lya luminosity function 
(LF) does not show any evolution between z — 3.1 - 5.7 
(Ouchi et al. 2008), it changes appreciably between z = 5.7 
and 6.6, with L, (the luminosity of the break, after which 
the number density decreases rapidly) at z = 6.6 (Kashikawa 
et al. 2006) being about 50% of the value at z — 5.7 (Shi- 
masaku et al. 2006). Unexpectedly, however, the ultraviolet 
(UV) LF does not show any evolution between these same 
redshifts. Although Kashikawa et al. (2006) have proposed 
this evolution in the Lya LF to be indicative of a sudden 
change in the ionization state of the IGM, the problem of 
why reionization would affect the high luminosity tail of the 
Lya LF rather than the faint end, as expected, remains. 

In spite of a number of diverse approaches, the effects 
of the ionization state of the the IGM on the visibility of 
LAEs still remain poorly understood. This is primarily be- 
cause understanding/ constraining the ionization state of the 
IGM using LAEs requires: (a) a detailed knowledge of the 
physical properties of each galaxy, including the SFR, stellar 
age and stellar metallicity, needed to calculate the intrinsic 
Lya and continuum luminosity produced by stellar sources, 
(b) the intrinsic Lya and continuum luminosity produced 
by the cooling of collisionally excited H i in the interstellar 
medium (ISM), (c) an understanding of the dust formation, 
enrichment and distribution in each galaxy, necessary to cal- 
culate the fractions of escaping Lya and continuum photons, 
and (d) a full radiative transfer (RT) calculation to obtain 
the fraction of Lya luminosity transmitted through the IGM 
for each galaxy. 

Although a number of simulations have been used in the 
past to study LAEs, they lack one or more of the aforemen- 
tioned ingredients. Since they used an N-body simulation, 
McQuinn et al. (2007) were forced to neglect the intrin- 
sic properties of the galaxies and their dust enrichment, al- 
though they carried out analytic RT. Iliev et al. (2008) used 
a simulation that only followed dark matter; they also could 
not use the intrinsic galaxy properties or calculate the dust 
enrichment, although they carried out a complete RT calcu- 
lation. Nagamine et al. (2008) used an SPH (smoothed par- 
ticle hydrodynamics) simulation to obtain the intrinsic SFR 
for each galaxy; however, they did not calculate the dust 
enrichment and assumed a fixed value of the IGM transmis- 
sion, ignoring RT. Dayal et al. (2009a, 2009b) used an SPH 
simulation and the intrinsic galaxy properties to calculate 
the luminosity from both stellar sources and cooling of H i , 
calculated the dust enrichment and the IGM transmission; 
the only missing ingredient in their work was the RT calcu- 
lation. Most recently, Zheng et al. (2009) have carried out 
a RT calculation on an SPH simulation; however, they have 
not included any dust calculation or the luminosity contri- 
bution from the cooling of H i . 

Our aim in this work is to build a LAE model con- 
taining all these ingredients, so as to determine the relative 
importance of dust, peculiar velocity fields and patchy reion- 
ization, in shaping the observed Lya and UV LFs. We start 
with an introduction to the SPH and RT simulations used 
in Sec. [5] We then present the main ingredients of our LAE 
model, including the calculation of the intrinsic luminos- 
ity from stellar sources/cooling of H i , the dust enrichment 



and the IGM transmission, in Sec. [3] Once the model is laid 
down, we present the results obtained with it and quan- 
tify the relative importance of dust, peculiar velocities and 
reionization on the Lya and UV LFs in Sec. [4] The physical 
properties of the galaxies identified as LAEs are shown in 
Sec. \S\ We conclude by mentioning the caveats and short- 
comings of our model in Sec. [6] The simulations used in 
this work are based on a ACDM cosmological model with 
fl A = 0.7, Q b = 0.3, n m = 0.04, H = WOh kms _1 Mpc _1 , 
where h — 0.7 and a scale invariant power spectrum of the 
initial density perturbations is normalized to as = 0.9. 



2 HYDRO AND RADIATIVE TRANSFER 

The Lya and UV LFs presented in this work are based 
on the results from combined runs of SPH and RT sim- 
ulations, carried out using GADGET-2 (Springel 20050 and 
CRASH (Maselli, Ferrara & Ciardi 2003; Maselli, Ciardi & 
Kanekar 2009) respectively, which are coupled to a previ- 
ously developed LAE model (Dayal et al. 2008; Dayal et 
al. 2009a, 2009b) . GADGET-2 generates the redshift evolution 
of the density field, the baryonic density distribution and 
the velocity fields in a lOO/i -1 Mpc comoving volume; we 
obtain a snapshot of the simulation at z ~ 6.1. The spe- 
cific simulation used in this work is the G5 run described 
in Springel and Hernquist (2003), which is part of an ac- 
curate study focused on modelling the star formation his- 
tory of the universe. We have chosen the G5 run since it 
contains several physical ingredients necessary for our in- 
vestigation. Firstly, the volume is large enough so that cos- 
mic variance is minimized in the determination of the LFs. 
Secondly, the SFR (iW*) is self-consistently inferred from 
physically motivated prescriptions that convert gas parti- 
cles into stellar particles and that properly take account of 
the mechanical/chemical feedback associated to supernovae 
and galactic winds. This point is particularly relevant since 
the intrinsic luminosity of the galaxies, both in the Lya and 
UV, depends sensitively on M» . Also, since the physics gov- 
erning star formation and galaxy evolution is modelled ac- 
curately, the simulation gives us a reliable representation 
of the galaxy population. As expected, the large simulation 
volume naturally leads to a relatively coarse mass resolution 
with the resolution mass being 2.12 x 10 10 Mq (3.26 x 10 s 
Mq) for dark matter (gas particles). Running a friends-of- 
friends (FOF) algorithm on the SPH particle distribution, 
we identify galaxies and obtain their intrinsic properties, 
including M t , the mass weighted stellar metallicity (Z»), 
the total gas and stellar masses. These are then used to 
calculate the intrinsic Lya/continuum luminosity, the dust 
enrichment and the escape fraction of continuum/Lya pho- 
tons from the galaxy. All these are presented and discussed 
in greater detail in Sec. [3] and [4] 

The RT calculations have been carried out in post- 
processing mode using the 3D CRASH code. The main as- 
sumptions made for running CRASH are the following: (a) 
we initially assume, at z ~ 6.1, the gas to be in photoion- 
ization equilibrium with a uniform ultraviolet background 
(UVB) produced by sources below the simulation resolution 

1 http:/ /www. mpa-garching.mpg.de/gadget/ 
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Figure 1. Maps (100/i — 1 Mpc on a side) of the spatial distribution of H I in a 2D cut through the RT simulation box showing the 
time evolution of the reionization process, with the colorbar showing the values (in log scale) of the H I plotted. The average decreasing 
neutral hydrogen fractions marked above each panel, (xm) = 0.295,0.24,0.157,4.45 X 10 _2 ,1.13 X 10~ 2 ,3.43 X 10 -3 , correspond to 
increasing RT simulation timescales of 10, 50, 100, 200, 300, 500 Myr, respectively. Details in text of Sec. [2] 



scale, corresponding to an average neutral hydrogen frac- 
tion, (xm) — 0.3, (b) we assume that galaxies have formed 
stars in the past at an average rate M*, equal to the one 
derived at z ~ 6.1 and (c) since the RT simulation is carried 
out on a SPH snapshot which contains about 2500 galaxies, 
to reduce the computational running time, we run the RT 
calculation in the monochromatic mode, with hv = 13.665 
eV. This is a reasonable assumption since the Lya LFs de- 
pend mostly on the average ionization fraction and on the 
3D topology of the fully ionized regions; we are therefore not 
interested in the detailed profiles of the ionization front or a 
highly accurate estimate of the IGM temperature, to probe 
which the ionizing radiation of each galaxy would have to 
be computed from its spectrum. 

The specific photoionization rate is calculated using the 
expression given in Eq. 12 of Dayal, Ferrara & Gallerani 
(2008) and then averaging this value over and stellar ages 
(t») of the galaxies in the the SPH simulation. The contri- 
bution from all the galaxies determined under the above- 
mentioned assumptions is then super-imposed on the av- 
erage UVB value. Finally, we use f esc = 0.02 (Gnedin et 
al. 2008), as a characteristic value of the escape fraction of 
H i ionizing photons from each galaxy. The total ionizing 
radiation field, which is the sum of that produced by all 
galaxies identified in the simulation and the UVB, is then 
evolved with CRASH up to the Hubble time corresponding 
to z ~ 6.1. However, as we are interested in assessing the 
visibility of LAE as a function of the mean ionization frac- 
tion which could correspond to different epochs along any 



assigned reionization history, no particular meaning should 
be attached to the redshift parameter. 

In Fig. [T]we show maps of the H i fraction across a 2D 
cut through the RT simulation box, at different simulation 
times. This time sequence then represents a time-line of the 
reionization process, due to the galaxy population present 
within the box. Though an approximation, this description 
catches many of the important features of the reionization 
process, including the complex topology produced by the 
inhomogeneities in the density field, by the galaxy properties 
and by their spatial distribution. 

Fig. [T]clearly shows the reionization process in its main 
three phases. First, isolated ionized hydrogen (H n ) regions 
start growing around point sources which are preferentially 
located along the over-dense filaments of matter. In this 
stage the size of the isolated H n regions grows differen- 
tially according to the production rate of H i ionizing pho- 
tons from the specific source; the few galaxies with large M* 
(e.g. the largest galaxy located at the center having M* ~ 
726 M0yr _1 ) can build an ionized region of about 20 Mpc, in 
a time as short as 10 Myr (see upper- left panel), which grows 
at a very fast rate due to the vigorous H i ionizing photon 
output. If observed at the same evolutionary timescale, due 
to the smaller H i ionizing photon budget, the H n regions 
of smaller galaxies are confined to sizes of the order of few 
Mpc (< 4 Mpc), which grow slowly. After 10 Myr of contin- 
uous star formation activity, the average neutral hydrogen 
fraction decreases only to (xhi) ~ 0.295. At this early stage, 
galaxies with M* ^ 25M0yr _1 , transmit only about 20% 
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of their Lya luminosity through the IGM. Consequently the 
low luminosity end of the LF is depressed with respect to 
the intrinsic emissivity of these sources. On the other hand, 
galaxies with larger M* are able to transmit a larger amount 
(~ 30%) of their luminosity through the IGM. 

As the H i ionizing photon production from galaxies 
continues, the sizes of the H n regions increase with time. At 
about 100 Myr (upper right panel), the isolated H n regions 
start overlapping, resulting in an enhancement of the lo- 
cal photo-ionization rate close to the smaller sources. This 
allows the ionized regions associated to the latter to grow 
faster, and results in the build-up of large H n regions also 
around smaller galaxies. This results in an enhancement in 
the transmission of the Lya luminosity for these sources, 
such that this value increases to ~ 40%, consequently lead- 
ing to a boost of the low-luminosity end of the Lya LF. As 
expected, the overlap phase speeds-up the reionization pro- 
cess; for a star formation time 50, 100, 200, 300 Myr, (xhi) 
decreases to ~ 0.24, 0.16, 4.5 x 10~ 2 , 1.1 x 10" 2 , respectively. 

Finally at about 300 Myr (lower center panel), the sim- 
ulation volume is almost completely ionized. This configu- 
ration corresponds to the post-overlap phase which leads to 
a further acceleration of the reionization process. For about 
400 Myr of star formation activity, (xhi) drops to a value 
of about 4.3 x 10 -3 . This is a consequence of the fact that, 
as the IGM becomes more ionized, even the ionizing ra- 
diation from the numerous smaller galaxies contributes to 
the overall radiation field. It is very important to point out 
that, even in this final stage, the ionization field is highly 
inhomogeneous, with regions of higher ionization fraction 
corresponding to the environment of point ionizing sources 
(green regions in the maps). This topology is of fundamen- 
tal importance in the estimation of the Lya LFs, as the Lya 
transmissivity associated to each source depends on the ex- 
tent of the highly ionized regions as well as on the residual 
neutral hydrogen within them; a very tiny fraction of resid- 
ual neutral hydrogen is sufficient to suppress the intrinsic 
Lya luminosity significantly. Even for an ionization fraction, 
(Xhi) ~ 3.4 x 10~ 3 , which results after 500 Myr of star for- 
mation, the transmission value ranges between 42 — 48% for 
galaxies with increasing M* ; the blue part of the line can be 
significantly attenuated even by the tiny residual fraction of 
about 10 -6 , close to the source. 



3 THE LAE MODEL 

As mentioned in Sec. [21 for each galaxy in the SPH simula- 
tion snapshot, we obtain the SFR, M„, the mass weighted 
stellar metallicity, Z* , the total gas mass, M g , and the stellar 
mass, M»; the stellar age is then calculated as t* = M*/M, , 
i.e. assuming a constant SFR. This assumption has been 
made for consistency with the RT calculation (see Sec. [2]). 
These properties are then used to calculate the intrinsic 
Lya (I/™') and continuum luminosity (L* nt ), considering the 
contribution from both stellar sources and from the cooling 
of collisionally excited H i in the ISM of the galaxy. The 
intrinsic values of the Lya and continuum luminosities so 
obtained can then be used to calculate the observed lumi- 
nosity values. It is comparatively easy to calculate the ob- 
served continuum luminosity, L c , since the wavelength range 
(1250-1500 A) of the continuum luminosity band is chosen 



so as to be unaffected by H i ; L c in fact only depends on 
the fraction of the continuum photons, f c , that escape the 
galactic environment undamped by dust, i.e. 

L c = 4 n Vc (1) 

On the other hand, since Lya photons have a large ab- 
sorption cross-section against H i , the observed Lya lumi- 
nosity, L a , depends both on: (a) the fraction of Lya photons, 
f a , that escape out of the galaxy into the IGM, undamped 
by dust and (b) the fraction, T a , of these "escaped" photons 
that are transmitted through the IGM, undamped by H i . 
The observed Lya luminosity can therefore be expressed as 

La = La faT a - (2) 

All galaxies with L a 10 42,2 ergs _1 and an observed 
equivalent width, EW = (L a /L c ) ^ 20 A are identified as 
LAEs and used to construct the cumulative Lya and UV 
LFs which can then be compared to the observations. In 
principle, our results can be compared to both the observed 
LFs at z ~ 5.7 and 6.6. However, there exists a huge un- 
certainty between the complete photometric and confirmed 
spectroscopic sample at z ~ 6.6; we limit the comparison 
between the theoretical model and observations to the data 
accumulated at z ~ 5.7 by Shimasaku et al. (2006), in this 
work. The calculations of L™', L z c nt , f a , f c and T a are ex- 
plained in more detail in what follows. 

3.1 The intrinsic luminosity model 

Star formation produces photons with energies > 1 Ryd, 
which ionize the H i in the ISM. Due to the high density of 
the ISM, recombinations take place on a short time scale and 
this produces a large number of Lya photons. The stellar lu- 
minosity component is calculated by using STARBURST99 
(Leitherer et al. 1999), a population synthesis code, to ob- 
tain the spectra and hence, the H i ionizing photon rate (Q) 
of each galaxy, using a Salpeter IMF, Z„, t* and M» as input 
parameters. Then the intrinsic Lya luminosity produced by 
stellar sources L a , is calculated as 
2 

L* a = ~(1 - fesc)Qhv a , (3) 

where f esc — 0.02, for consistency with the RT calculation 
(see Sec. [2}, h is the Planck constant, u a is the frequency 
of the Lya photons and the factor of two-thirds enters as- 
suming case B recombination (Osterbrock 1989) . The intrin- 
sic continuum luminosity produced by stellar sources, L* , is 
calculated at 1375 A, at the middle of the continuum band 
between 1250-1500 A. 

The Lya and continuum luminosity from the cooling 
of collisionally excited H i in each galaxy (L a ,L^. respec- 
tively) is calculated using the number density of electrons 
and H i , which depend on the temperature distribution 
of the ISM gas mass and the gas distribution scale of the 
galaxy. The temperature dependence comes from the fact 
that the recombination coefficient, which determines the 
neutral hydrogen fraction in the ISM, decreases sensitively 
with higher temperatures. The H i number density is cal- 
culated by assuming a primordial ISM composition (76% 
hydrogen and 24% helium) and calculating the gas distribu- 
tion scale assuming the H i to be concentrated in a radius, 
r g = 4.5Ar2oo- Here the spin parameter, A — 0.05 (Ferrara, 
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Pettini & Shchekinov 2000) and V2oo is the virial radius, as- 
suming the collapsed region has an overdensity of 200 times 
the critical density at the redshift considered. The calcula- 
tion of the electron number density requires a knowledge of 
the temperature distribution of the IGM gas. However, since 
we do not have this information from the SPH simulation 
used in this work, we use the temperature distribution av- 
eraged over distinct halo mass ranges shown in Tab. 1 of 
Dayal et al. (2009b) and again use the gas distribution scale 
to calculate the electron number density. Complete details 
of this calculation can be found in Sec. 3.2 of Dayal et al. 
(2009b). 



3.2 The dust model 

It is important to model the dust in the ISM of galaxies 
since it absorbs both Lya and continuum photons, thereby, 
strongly influencing the observed luminosity values. Dust is 
produced both by supernovae and evolved stars in a galaxy. 
However, several authors (Todini & Ferrara 2001; Dwek et 
al. 2007) have shown that the contribution of AGB stars can 
be neglected for z ^ 5.7, since the typical evolutionary time- 
scale of these stars (^ 1 Gyr) becomes longer than the age 
of the Universe above that redshift. Hence, assuming SNII 
to be the primary dust factories at z ~ 6.1, we calculate the 
dust enrichment for each galaxy, taking into account three 
processes: (a) SNII produce dust in the expanding ejecta; 
the average dust mass produced per SNII is taken to be 
0.54 M (Todini & Ferrara 2001; Nozawa et al. 2003, 2007; 
Bianchi & Schneider 2007), (b) SNII destroy dust in the 
ISM they shock to velocities ^ 100 km s -1 , for which we 
use a destruction efficiency value ~ 0.12 (McKee 1989), and 
(c) a homogeneous mixture of gas and dust is assimilated 
into star formation (astration). Once the final dust mass, 
Mdust{t*), is calculated for each galaxy in the simulation, we 
can translate this into an optical depth, r c , for continuum 
photons as 



3£d 



(4) 



where Ed = Md us t{t»)r^ 2 is the dust surface mass density, 
rd is the radius of dust distribution, a and s are the ra- 
dius and material density of graphite/carbonaceous grains, 
respectively (a = 0.05/im, s = 2.25 gem -3 ; Todini & Fer- 
rara 2001; Nozawa et al. 2003). This can be used to obtain 
the escape fraction of continuum photons from the galaxy 
as f c — e -T< = . The determination of rd follows in Sec. 0] and 
complete details of the calculations mentioned here can be 
found in Dayal et al. (2009b). 

Different extinction curves, including that for the Milky 
Way, Small Magellanic Cloud and for supernovae give differ- 
ent relations between the extinction of the Lya and contin- 
uum photons, for homogeneously distributed dust. In this 
spirit, and to get a hint of the dust inhomogeneity, we com- 
bine f a and f c , and use the escape fraction of Lya pho- 
tons relative to the continuum ones, fa./ fc, to be the free 
parameter to calculate L a in our model. A more detailed 
explanation of the determination of f c and f a , required to 
reproduce the observations, follows later in Sec. [4] 



3.3 The IGM transmission model 

The calculation of T a is now explained in better detail. If 
z em an d Zobs are the redshifts of the emitter and the observer 
respectively, we calculate the total optical depth (r a ) to the 
Lya photons along a line of sight (LOS) as 



T a (v) 



<Ta<t>(v)nHi(z)^dz, 



(5) 



where v = (A — \ a )[\ a c\~ 1 is the rest-frame velocity of a 
photon with wavelength A, relative to the line centre (rest- 
frame wavelength X a = 1216 A, velocity v a = 0), (f> is the 
Voigt profile, uhi is the number density of H i along the 
LOS, dl/dz = c[H(z)(l + z)]' 1 and c represents the light 
speed. We use ao = 7re 2 /[m e c] _1 , where e, m e are the elec- 
tron charge and mass respectively and / is the oscillator 
strength (0.4162). 

For regions of low H i density, the natural line broad- 
ening is not very important and the Voigt profile can be 
approximated by the Gaussian core: 



Xa 



-\ %- 



-r 2 



(6) 



where v; is the velocity of a photon of initial velocity v at 
a redshift z; along the LOS, v p is the peculiar velocity at Zi 
and v a is the velocity of the Lya photons at z%, which is in 
this expression. The Doppler width parameter is expressed 
as b = y^2kT/m,H , where mu is the hydrogen mass, k is the 
Boltzmann constant and T is the IGM temperature. While 
we use a T = 10 4 K for each cell with an ionization fraction 
(1 - Xhi) > 0.1, we use T = 20 K if the IGM is pristine, i.e. 
has an ionization fraction (1 — xhi) ^ 0.1. Q 

In more dense regions the Lorentzian damping wing of 
the Voigt profile becomes important. According to Peebles 
(1993), this can be approximated as 



^lorentz (v?) 



R,aX Q 



7T[(vi + v p — v a )' 2 + R a ] 



(7) 



where R a = AA Q [47r] _1 and A = 6.25 x 10 s s _1 is the decay 
constant for the Lya resonance. Although computationally 
more expensive than the above approximations, using the 
Voigt profile to compute the absorption cross-section gives 
precise results, and therefore we have implemented it in our 
code to obtain all the results presented below. 

The total Lya optical depth is calculated for each 
galaxy identified in the CRASH outputs corresponding to 
(Xhi) ~ 0.29,0.24,0.16,4.5 x 10 _2 ,1.1 x 10 _2 ,4.3 x 
10 _3 ,3.4 x 10 -3 as mentioned in Sec. [2] The transmission 
along any LOS is calculated as T a = e~ T ° , integrating from 
the position of each galaxy to the edge of the box. To get a 
statistical estimate of the transmission, we construct 6 LOS 
starting from the position of each galaxy, one along each 
of the box axes, for each of the mentioned ionization states. 
The average value of the transmission for each galaxy is then 
obtained by averaging the transmission over the 6 LOS. This 



2 We find that even for an average neutral hydrogen fraction 
as high as, (xhi) ~ 0.295, only 8 cells out of 256 3 have (1 — 
Xhi) ^ 0.1, with the complete RT simulation as described in 
Sec. [2] Hence, to simplify the calculations, we use T = 10 4 K for 
all the cells for all the CRASH outputs. 
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Figure 2. Cumulative Lye? LF for galaxies identified as LAEs 
from the SPH simulation snapshot including the full RT calcula- 
tion and density fields but ignoring velocity fields (v p = 0) and 
dust {f a = / c = 1). The lines from bottom to top correspond to 
CRASH outputs with decreasing values of (xhi) such that: (xhi) ~ 
= 0.29 (solid), 0.24 (dotted), 0.16 (short dashed), 4.5 X 10~ 2 
(long dashed), 1.1 X 10~ 2 (dot-short dashed), 4.3 X 10~ 3 (dot- 
long dashed) and 3.4 X 10 -3 (short-long dashed). Points show 
the observed LF at z ~ 5.7 (Shimasaku et al. 2006). 



is done for all the ionization state configurations resulting 
from the RT calculation. 

Once the above model is in place, the only two free pa- 
rameters of our model to match the theoretical Lya and 
UV LFs to the observed ones are: the dust distribution ra- 
dius, rd, and the relative escape fraction of Lycv photons 
as compared to the continuum photons, f a /fc- Both these 
parameters remain quite poorly understood due to a lack 
of observational data about the dust distribution/topology 
in high-redshift galaxies; they must therefore be inferred by 
comparing the theoretical LFs to the observed ones. Details 
on how these parameters are determined follow in Sec. [4] 



4 LAE VISIBILITY DURING REIONIZATION 

Once the combined SPH+RT calculations are carried out 
and the LAE model implemented, we are in a position to 
quantify the importance of reionization, peculiar velocities 
and the dust enrichment on the observed LFs. Physically, 
the reionization process leads to a decrease in (xhi), thereby 
increasing T a ; on the other hand, peculiar velocities caused 
by galactic scale outflows (inflows), redshift (blueshift) the 
Lycv photons, thereby leading to a higher (lower) value of 
T a . A handle on the dust enrichment is necessary since dust 
grains absorb both Lya and continuum photons, thereby 
affecting their escape fractions from the galaxy. 

We begin our study by ignoring the effects of the pecu- 
liar velocities and assuming all the galaxies to be dust free, 
to quantify how each of these two parameters shapes the 
observed Lya and UV LFs. The former assumption means 
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Figure 3. Cumulative Lye? LF for galaxies identified as LAEs 
from the SPH simulation snapshot including the full RT calcula- 
tion, density fields and dust but ignoring velocity fields (v p = 0), 
corresponding to models SI — S7 in Tab. [TJ The lines from bot- 
tom to top correspond to CRASH outputs with decreasing values of 
(XHl) such that: (xsi) ~ 0.29 (solid), 0.24 (dotted), 0.16 (short 
dashed), 4.5 X 10 -2 (long dashed), 1.1 X 10 -2 (dot-short dashed), 
4.3 X 10~ 3 (dot-long dashed) and 3.4 X 10 -3 (short-long dashed). 
The shaded region shows the poissonian error corresponding to 
(Xhi) ~ 3.4 X 10 -3 , i.e., model S7 (Tab.[TJ and points show the 
observed LF at 2 ~ 5.7 (Shimasaku et al. 2006). 



that v p = in Eqs. [6] [7] the latter implies that all pho- 
tons produced inside the galaxy escape into the IGM, i.e. 
fc = fa = 1 in Eqs. [1] [2] We then use the prescriptions de- 
tailed in Sec.[3]to identify the galaxies that would be visible 
as LAEs in each CRASH output and build their Lya LF, as 
shown in Fig. [2] 

Since both peculiar velocities and dust are neglected, 
in this case, while the UV LF is simply the intrinsic LF, 
the Lya LF is shaped solely by the transmission through 
the IGM. As mentioned in Sec. [2l after a star formation 
time scale as short as 10 Myr ((xhi) ~ 0.29), the galax- 
ies in the simulation snapshot are able to form H n regions, 
with the region size and the ionization fraction inside it in- 
creasing with M»; this results in a T a which increases with 
M*. As the galaxies continue to form stars, the sizes of the 
H ii regions built by each source increase with time, leading 
(Xhi) to decrease to 0.24,0.16 at t = 50,100 Myr respec- 
tively. This leads to an increase in the transmission of the 
red part of the Lya line for all the sources, in turn yielding a 
corresponding increase of the Lya luminosity of each source, 
as depicted in Fig. [2] However, after about 200 Myr, (xhi) 
reduces to ~ 0.04; the H n regions built by each source are 
large enough so that almost all of the red part of the Lya line 
is transmitted and the value of T a saturates for all LAEs. 
This results in very similar LFs for (xhi) ^ 0.04, as seen 
from the same figure. 

The key point here is that, if peculiar velocities and dust 
effects are neglected, at no stage of the reionization history, 
either the slope or the amplitude of the observed Lya LF 
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Figure 4. T a as a function of M, for galaxies identified as LAEs from the simulation snapshot including the full RT calculation, 
density fields, dust and ignoring (including) velocity fields shown by dotted (solid) lines in each panel. The solid/dotted lines in each 
panel correspond to the following models given in Tab. [TJ (a) Ml/Sl {{xhi) ~ 0.29), (b) M2/S2 (0.24), (c) M3/S3 (0.16), (d) M4/S4 
(4.5 X 10~ 2 ), (e) M5/S5 (1.1 X 10~ 2 ) and (f) M7/S7 (3.4 X 10~ 3 ). In each panel, the stellar mass bins span 0.2 dcx and the shaded 
regions represent the la error bars in each mass bin. 



can be reproduced; analogous problems arise also when the 
UV LFs are considered. 

The above discussion implies the need of one or more 
physical effects attenuating the Lya and continuum photons, 
the most obvious of which is the presence of dust in the 
ISM of these galaxies, which would absorb both Lya and 
continuum photons, simultaneously reducing f a and f c . We 
then include the dust model described in Sec. 13.21 into our 
calculations. Since the UV is unaffected by the H i in the 
IGM as mentioned in Sec. [3] the same value of r& — 0.48r 9 
(r g is the ISM gas distribution scale, Sec. l3.2|) reproduces the 
UV LF for all of the ionization states of the IGM, ranging 
from {xhi) ~ 0.29 to 3.4 X 10" 3 . This assumption of the 
dust distribution scale leads to an average escape fraction of 
continuum photons, f c ~ 0.12 (refer Sec. 13.21 for details) for 
the galaxies we identify as LAEs. 

We now momentarily digress to discuss the effect of 
reionization on T a (defined by Eq. [2| for different stellar 
mass ranges. However, since M* scales with M» (details in 
Sec. [5]), we discuss T a in terms of M*. As mentioned in Sec.[2l 
we initialize the CRASH runs with (xhi) ~ 0.3. In a star for- 
mation timescale of 10 Myr, by virtue of the H n regions 
that already start growing, the total photoionization rate 
(sum of contributions from the UVB and all the galaxies 
in the simulation) has a value of Ft0 = 2.8 X 10 _17 s _1 
and (xhi) decreases slightly to 0.295. At this point, due 



3 Tt, the total photoionization rate, is calculated assuming 
ionization-recombination equilibrium over average values of den- 
sity and ionization fraction in the simulation volume, so that 

r = (1 ~ {XHl)) 2 {n H )a B 
{Xhi) 



to their smaller H i ionizing photon output, galaxies with 
M* ^ 25Mgyr , have T a ~ 0.2; galaxies with larger M*, 
have T a ~ 0.3, as shown in Panel (a) of Fig. [4] As the star 
formation continues, for 50 (100) Myr of star formation, 
the H ii region sizes increase and Ft increases slightly to 
~ 4.0 x 10~ 17 (7.7 x 10~ 17 ) s" 1 ; T a increases and ranges be- 
tween 0.28-0.4 (0.34-0.44) for M» ~ 8-200 M yr" 1 , Panel 
b (c), Fig.Q] Finally for Ft ~ 3.3 x 10" 16 s _1 , corresponding 
to {xhi) ~ 0.04, the transmission settles to T a ~ 0.4 — 0.48 
for M* ~ 8- 200M©yr _1 (Panels d-f); in about 200 Myr 
from the ignition of star formation, the Stromgren spheres 
built by these LAEs are large enough so that the redshifted 
Lya photons are no longer affected by the H i outside this 
region. However, the residual H i inside this ionized region 
leads to an absorption of the photons blueward of the Lya 
line, and hence, about half of the line is transmitted. The 
values of Ft, {xhi) and the T a values averaged for all LAEs 
in each CRASH output are shown in Tab. [T] 

Using the above T a values and dust model fixed by the 
UV LF, the only free parameter we are left with, to match 
the theoretical and observed Lya LFs is f a /fc', this only 
scales the Lya LF without affecting its shape. As mentioned 
above, the T a value for galaxies identified as LAEs settles 
for {xhi) ^ 0.04, which leads to a corresponding saturation 
in the Lya LF. We find that f a /fc ~ 1.3 reproduces the 
slope and the magnitude of the observed Lya LF quite well 
for all the CRASH outputs where {xhi) ^ 0.04, as shown in 
Fig. |31 However, due to decreasing values of T a , the Lya 



Here, (njj) is the average hydrogen number density in the simu- 
lation volume and as is the case B recombination co-efficient. 
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Figure 5. UV LF for galaxies identified as LAEs from the simula- 
tion snapshot including the full RT calculation, density fields, ve- 
locity fields and dust, corresponding to models Ml to M7 in Tab. 
[T] The lines from bottom to top correspond to decreasing values of 
(Xhi) such that: (xhi) ~ 0.29 (solid), 0.24 (dotted), 0.16 (short 
dashed), 4.5 X 10~ 2 (long dashed), 1.1 X 10~ 2 (dot-short dashed), 
4.3 X 10~ 3 (dot-long dashed) and 3.4 X 10 -3 (short-long dashed). 
The shaded region shows the poissonian error corresponding to 
(Xhi) ~ 3.4 X 10~ 3 , i.e., model M7 (Tab [TJ and points show the 
observed UV LF at z ~ 5.7 (Shimasaku et al. 2006). 
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Figure 6. Cumulative Lya LF for galaxies identified as LAEs 
from the simulation snapshot including the full RT calculation, 
density fields, velocity fields and dust, corresponding to models 
Ml to M7 in Tab. [Tl The lines from bottom to top correspond to 
decreasing values of (xni) such that: (xhi) ~ 0.29 (solid), 0.24 
(dotted), 0.16 (short dashed), 4.5 xlO -2 (long dashed), l.lxlO -2 
(dot-short dashed), 4.3 X 10 -3 (dot-long dashed) and 3.4 X 10 -3 
(short- long dashed). The shaded region shows the poissonian error 
corresponding to (xhi) ~ 3.4 X 10~ 3 , i.e. model M7 (Tab.QJ and 
points show the observed LF at 2 ~ 5.7 (Shimasaku et al. 2006). 



LF is progressively under-estimated for increasing values of 
{Xhi}, as seen from the same figure. 

Once the above framework is in place, we include and 
study the effect of the final component that can affect the ob- 
served luminosity, the presence of peculiar velocities. Once 
this is done by consistently deriving the peculiar velocity- 
field from the simulations, we again reproduce both the mag- 
nitude and slope of the observed UV LF as shown in Fig. [5] 
Of course this does not come as a surprise since velocity 
fields do not affect UV photons. However, galactic scale out- 
flows (inflows) from a galaxy, redshift (blueshift) the Lya 
photons, thereby leading to a higher (lower) T a value. It is 
quite interesting to see that we can again match the the- 
oretical Lya LFs to the observed ones by a simple scaling 
between f a and f c for {xhi} ^ 0.04, as seen from Fig. [5] 
However, we require a much higher value of f a /fc ~ 3.7 
(models M4 — M7, Tab.[l}, compared to the value of 1.3 ex- 
cluding velocity fields (models 54 — S7, Tab. [[J. Again, T a , 
and hence the Lya LF, get progressively more damped with 
an increase in {xhi}, as seen from Fig. [6] (models M3 — Ml, 
Tab. [J). 

We now discuss the reason for the higher fa/fc ~ 3.7 
value required to reproduce the Lya LF when velocity fields 
are considered, as compared to the ratio of 1.3, when they 
are not, for (xhi) ^ 0.04. The galaxies we identify as LAEs 
have halo masses Mh ~ igio.4-12 jy^ Q] which correspond to 
^ 2cr fluctuations at z ~ 6.1. These LAEs are therefore sub- 
ject to strong inflows since they lie in dense regions. As these 
inflows blue-shift the Lya photons, even photons in the red 
part of the line are attenuated, thereby reducing T a . This 



can be seen clearly from Fig.[4l where the solid lines, which 
represent the outcomes from the models including velocities, 
are always significantly below the dotted ones, which repre- 
sent the model with no velocity field included. When velocity 
fields are included we find T a ~ 0.08 - 0.12 (0.16 - 0.18) for 
{Xhi) ~ 0.29 (3.4 X 10~ 3 ) for M» ~ 8 - 200 M yr" 1 . Aver- 
aging over all the LAEs for {xhi) ~ 3.4 X 10 -3 , we find the 
value of T a ~ 0.17 (0.44) including (excluding) the effects 
of peculiar velocity fields. This requires that when pecu- 
liar velocities are included, a correspondingly larger fraction 
(0.44/0.17 = 2.6) of Lya photons must escape the galaxy, 
undamped by dust to bring the observed luminosity up to 
the levels it would reach in the absence of these inflows. Al- 
though many galaxies do show outflows, powered by super- 
nova explosions, these dominate at very small scales (^ 170 
physical Mpc). However, the small redshift boost imparted 
by these is negligible compared to the blue-shifting of the 
Lya line because of large scale inflows; eventually, it is these 
dominant inflows that determine the value of T a . 

Further, including velocity fields changes the slope of 
the LF; when velocity fields are not included, T a basically 
scales with M„; when these are included, T a is the most 
damped for the largest masses, since these see the strongest 
inflow velocities by virtue of their largest potential wells. 
This can be seen clearly from Fig. [4] where the slope of T a 
is visibly shallower for the solid curves (including velocities) 
than for the dotted ones (neglecting velocities) . This has the 
effect of flattening the slope of the Lya LF as can be seen 
in Fig. [6)1. 

As an important result, the above analysis shows that 
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Table 1. The model designation (col 1), the star formation timescale to ionize the IGM (col 2), the features included in the model (VF 
stands for peculiar velocity fields) (col 3), the average neutral hydrogen fraction (col 4), the value of the total photoionization rate (sum 
of the contribution from the UVB and from all the galaxies in the snapshot) corresponding to this neutral hydrogen fraction (col 5), the 
average transmission for all the galaxies identified as LAEs (col 6), the average escape fraction of continuum photons for LAEs for the 
model presented in Sec. 13.21 (col 7), the average color excess value corresponding to this continuum escape fraction using the supernova 
extinction curve (col 8) and the relative escape fraction of Lyc? and continuum photons to best fit the observations for (xhi) ^ 0.04 (col 
9). 
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Figure 7. The 1 — 5<r probability contours (black to light gray 
respectively) for combinations of {xhi) and fa/fc that fit the 
observed Lya LF. Within a la (5<r) error, we can not distinguish 
an IGM with (xhi) ~ 3.4 X 10 — 3 from one where (xhi) ~ 0.16 
(0.24) for faff c ~ 3.4 - 4.1 (3.0 - 5.7). This implies a degener- 
acy between the ionization state of the IGM and dust clumping 
(or grain properties) inside high-redshift galaxies; the ionization 
state of the IGM cannot be tightly constrained unless the rela- 
tive escape fraction of Lyo compared to the continuum photons 
is reasonably well understood. Refer Sec. [4] for details. 



there exists a degeneracy between the ionization state of the 
IGM and dust clumping (or grain properties) inside high- 
redshift galaxies, i.e a high (low) T a can be compensated 
by a low (high) f a . This is shown in Fig. [7] where we find 
that within a la error, for f a /f c ~ 3.4 — 4.1, we can not 
distinguish an IGM with (xhi) ~ 3.4 x 10 -3 from one where 



(Xhi) ~ 0.16. Within the area under the 5a error, fa/fc ~ 
4.1-5.7 can also fit the Lya LF for (xhi) ~ 0.24. This leads 
to the very interesting conclusion that the ionization state 
of the IGM cannot be tightly constrained unless the relative 
escape fraction of Lya compared to the continuum photons 
is reasonably well understood. 

We now discuss the dusty nature of the LAEs we iden- 
tify from the SPH simulation. We require that for the com- 
plete LAE model (Ml - Ml, Tab. the value of f a /f c 
must range between 3.4 — 4.1 (3 — 5.7) for an average neu- 
tral hydrogen fraction of (xhi) ^ 0.16 0.24). However, 
no single extinction curve gives a value of f a /fc > 1. One 
of the simplest ways of explaining this large relative escape 
fraction is to invoke the multiphase ISM model as proposed 
by Neufeld (1991), wherein the ISM is multiphase and con- 
sists of a warm gas with cold dust clumps embedded in it. 
This inhomogeneity of the dust distribution can then lead 
to a larger attenuation of the continuum photons relative to 
the Lya. 

We also translate f c into the color excess, E(B — V), 
for models Ml — M7 (Tab. [1} and compare this value to 
other high redshift LAE observations. At high redshifts, the 
observed properties of the most distant quasars (Maiolino 
et al. 2006) and gamma-ray bursts (Stratta et al. 2007) can 
be successfully interpreted using a SN extinction curve (To- 
dini & Ferrara 2001; Bianchi & Schneider 2007). Using the 
same curve, we find the average value of the color excess, 
E(B — V) = 0.2, corresponding to an average f c = 0.12. 
This inferred color excess value is in good agreement with 
recent experimental determinations: by fitting the SEDs of 3 
LAEs at z = 5.7, Lai et al. (2007) have inferred E(B-V) < 
0.225 — 0.425; in a sample of 12 LAEs at z = 4.5, Finkelstein 
et al. (2009) have found E(B - V) = 0.035 - 0.316; finally, 
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Figure 8. Physical properties of the galaxies identified as LAEs using model M7 in Tab.[TJ This model includes the full RT calculation, 
density, velocity fields, dust and has (xhi) ~ 3.4 X 10 -3 . The panels show (a) the gas mass, M g (b) SFR, M*, (c) mass weighted stellar 
age, t t and (d) the mass weighted stellar metallicity, Z t , (e) the total dust mass, M dust and, (f) the escape fraction of Lyc? photons, f a , 
plotted as a function of the stellar mass, M*. The values of each quantity are shown averaged over M* bins spanning 0.2 dex, with the 
error bars representing the la error in each bin. 



Pirzkal et al. (2007) have found E{B - V) = 0.025 - 0.3 for 
3 galaxies at z — 4 — 5.76. 

We now summarize our two main results: (a) we find 
that the Lya LF can be well reproduced (to within a 5a 
error) by an average neutral hydrogen fraction as high as 
0.24 (an almost neutral IGM), to a value as low as 3.4 x 
10 -3 , corresponding to an ionized IGM, provided that the 
increase in the transmission is compensated by a decrease in 
the Lya escape fraction from the galaxy, (b) we find that to 
reproduce the Lya LF, for any ionization state of the IGM, 
we require fn/fc > 1, a value that cannot be obtained using 
any existing extinction curve; this raises the need to invoke a 
multiphase ISM model, in which dust clumps are embedded 
in a highly ionized ISM, to facilitate the Lya photon escape 
relative to that of continuum photons. 



5 PHYSICAL PROPERTIES OF LAES 

We are now in a position to discuss the physical properties 
of the galaxies we identify as LAEs, including the stellar 
and gas mass, metallicity, stellar ages, dust mass and es- 
cape fraction of Lya photons from the ISM. As mentioned 
in Sec. [4] while the UV LF is independent of the ionization 
state of the IGM, the Lya LF can be well reproduced if a low 
(high) Lya escape fraction from the galaxy is compensated 
by a high (low) transmission through the IGM. Hence, by 
scaling up f a /fc for increasing values of (xhi) > 0.04 (mod- 
els M1-M3, Tab. [TJ, we would broadly always identify the 
same galaxy population as LAEs. We now show the physi- 
cal properties for the LAEs identified using model M7 (Tab. 



[TJ, which includes the full RT calculation, density /velocity 
fields and dust, with f a / fc = 3.7. 

We find that there is the direct correlation between M g 
and M*; LAEs with a larger M« are also more gas rich. 
However, the ratio M g /M* is the largest (smallest) for LAEs 
with the smallest (largest) M* , as one can deduce from panel 
(a) of Fig. [8] If we assume the halo mass to scale with the 
total baryonic mass (M g + M*), according to the cosmolog- 
ical ratio, this implies that more massive galaxies are more 
efficient in turning their gas into stars. This trend is as ex- 
pected, since even a small amount of star formation activity 
in low mass galaxies can lead to large outflows of gas, thereby 
suppressing further star formation. This situation however, 
does not occur in galaxies with larger masses, which do not 
witness large, galactic scale outflows, by virtue of their much 
larger potential wells (Mac Low & Ferrara 1999). 

LAEs with Af* < 10 9 - 7 M© have M* ~ 8 - 25 M yr~\ 
i.e. most of the LAEs have a sustained but not exceptionally 
large star formation activity. However, galaxies with larger 
stellar masses are much more efficient in converting their gas 
into stars, leading to SFR as large as 200 Mq yr _1 , as seen 
from panel (b), Fig. [5] which is just a positive feedback of 
the M g — M* relation mentioned above. 

Although the average age of the galaxies are calculated 
as t, = M*/M«, it is interesting to see that all the galaxies 
we identify as LAEs have t* ^ 10 Myr, as shown in panel (c) 
of Fig. [8] further the standard deviation on the smallest ages 
are the smallest. Since we calculate the ages assuming M» to 
have been constant over the entire star formation history, it 
is entirely possible that many of the LAEs are actually older 
(younger) if M, in the past was smaller (larger) than the 
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final value at z ~ 6.1. However, we are unable to comment 
on this further in absence of the complete star formation 
history for each galaxy. Using these average ages, we find 
that LAEs are intermediate age galaxies, instead of being 
extremely young (< 10 Myr) or extremely old (~ 1 Gyr) 
objects. This confirms an earlier result obtained in Dayal 
et al. (2009a), based on an accurate modelling of the star 
formation history, that it is unlikely that young ages could 
be responsible for the large equivalent widths (EW ^ 200 A) 
observed for a number of LAEs at z ~ 4.5, 5.7 by Dawson 
et al. (2007) and Shimasaku et al. (2006) respectively. A 
complete discussion on the EW distribution is deferred to 
Sec. [6] of this paper. 

The mass weighted stellar metallicity of LAEs scales 
with M* as shown in panel (d) of Fig. [8] LAEs with a higher 
M» are more dust enriched, which is only to be expected 
since the metals have a stellar origin. The metallicity values 
for the smallest halos show the largest dispersion, possibly 
arising due to the differing values of feedback in these low 
mass halos. Compared to the analogous mass-metallicity re- 
lation observed at lower redshifts (Tremonti et al. 2004; Pan- 
ter et al. 2008; Maiolino et al. 2008), we do not see the sign 
of a flattening of metallicity towards larger masses which 
might imply that at high redshift galaxies are not massive 
enough to retain their metals and behave as close-boxes. 

In our model, we have assumed SNII to be the primary 
dust factories, with the SNII rat^B being proportional to 
M*. In our model, the dust amount is regulated solely by 
stellar processes: dust is produced by SNII, destroyed in the 
ISM shocked by SNII and astrated into stars, as mentioned 
before in Sec. 13.21 The dust amount, therefore, scales with 
M* , with the most star forming galaxies being most dust en- 
riched, as shown in panel (e) of Fig.[H] A caveat that must be 
mentioned here is that as outflows tend to occur on smaller 
scales with respect to inflows, they are only marginally re- 
solved by our RT simulations. Lacking this information, we 
have not included their destructive impact into the compu- 
tation of the dust mass, which could therefore be somewhat 
overestimated (refer Dayal et al. 2009b for details). 

As expected, we find that f a decreases with increasing 
dust enrichment of the galaxy; galaxies with larger M* , and 
hence, M„, have a smaller Lya escape fraction. The value of 
f a decreases by a factor of about 3, going from 0.9 to 0.3 as 
M* runs from io 8 - 5 " 10 - 7 M©, as shown in panel (f), Fig. [SJ 



6 CONCLUSIONS AND DISCUSSION 

Although a number of simulations have been used in the 
past few years (Mc Quinn et al. 2007; Iliev et al. 2008, 
Nagamine et al. 2008; Dayal et al. 2009a, 2009b; Zheng et 
al. 2009), the relative importance of reionization, velocity 
fields and dust, on shaping the observed Lya and UV LFs 
has remained rather obscure so far. To this end, we build a 
coherent model for LAEs, which includes: (a) cosmological 
SPH simulations run using GADGET-2, to obtain the galaxy 
properties (M„, Z, , M*), (b) the Lya and continuum lumi- 
nosities produced, calculated using the intrinsic properties of 

4 The SNII rate is estimated to be 7M*, where 7 ~ (54 M ) -1 
for a Salpeter IMF between the lower and upper mass limits of 1 
and 100 Mq respectively (Ferrara, Pettini and Shchekinov 2000). 



each galaxy, accounting for both the contribution from stel- 
lar sources and cooling of collisionally excited H 1 in their 
ISM, (c) a model for the dust enrichment, influencing the 
escape fraction of Lya and continuum photons, and (d) a 
complete RT calculation, carried out using CRASH, including 
the effects of density and velocity fields, used to calculate 
the transmission of Lya photons through the IGM. 

In spite of this wealth of physical effects, we are still 
left with two free parameters: the radius of dust distribu- 
tion, Td, which determines the optical depth to continuum 
photons (see Eq. [4| and, the escape fraction of Lya photons 
relative to the continuum photons, f a / f c . These then must 
be constrained by comparing the theoretical model to the 
observations. 

Starting our analysis by assuming all galaxies to be 
dust-free and ignoring gas peculiar velocities, we find that 
with such a scheme, none of the ionization states of the 
IGM with (xhi) ~ 0.29 to 3.4 x 10~ 3 , can reproduce either 
the slope or the galaxy number density of the observed UV 
or Lya LF. We then include the dust model described in 
Sec. 13.21 into our calculations, such that each SNII produces 
0.54 Mq of dust, SNII shocks destroy dust with an efficiency 
of about 12% and the dust distribution radius is about half 
of the gas distribution radius, ra — 0.48r 9 . Since the UV is 
unaffected by the H 1 in the IGM as mentioned in Sec. [3] the 
same dust model parameters reproduce the UV LF for all of 
the ionization states of the IGM, ranging from (\hi) ~ 0.29 
to 3.4 x 10" 3 . 

We find that we can reproduce the Lya LF using 
fa/ fc ~ 1-3 (3.7) for all ionization states such that (xhi) ^ 
0.04 excluding (including) peculiar velocity fields, since the 
Lya LF settles to a constant value at this point due to a 
saturation in T a , as explained in Sec. 0] Higher values of 
(Xhi) naturally lead to a lower transmission, thus leading 
to a progressive underestimation of the Lya LF. The higher 
fa/ fc value required to reproduce the observations when ve- 
locity fields are included, arise because LAEs reside in ^ 2a 
fluctuations; LAEs are therefore subject to strong inflows 
which blue-shift the Lya photons, thereby reducing T a . This 
decrease in T a must therefore be compensated by a larger 
escape fraction from the galaxy to fit the observations. Fur- 
ther, due to the larger masses, T a is most damped for the 
largest galaxies, when velocity fields are included, which has 
the effect of flattening the slope of the Lya LF. 

The above degeneracy between the ionization state of 
the IGM and the dust distribution/clumping inside high- 
redshift galaxies has been quantified (see Fig. [TJl ; the Lya 
LF can be well reproduced (to within a 5a error) by (xhi) ~ 
0.24, corresponding to a highly neutral IGM, to a value as 
low as 3.4 x 10 -3 , corresponding to an ionized IGM, pro- 
vided that the increase in T a is compensated by a decrease 
in the Lya escape fraction from the galaxy. This leads to 
the very interesting conclusion that the ionization state of 
the IGM can not be constrained unless the escape fraction 
of Lya versus continuum photons is reasonably well under- 
stood. This is possible only through simulations that model 
the physics of the ISM, including the clumping of dust, and 
include a RT calculation for both Lya and continuum pho- 
tons on such ISM topologies. However, the clumping of dust 
depends on the turbulence in the ISM, which itself remains 
only poorly understood. Therefore, attempts to break the 
mentioned degeneracy must be deferred to future works. 
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As for the dusty nature of LAEs, we find that (f c ) ~ 
0.12, averaged over all the LAEs in any snapshot, which 
corresponds to a color excess, E(B — V) ~ 0.2, using a SN 
extinction curve. This value is in very good accordance with 
the observed values of E(B -V) < 0.225 - 0.425 (Lai et al. 
2007), E(B - V) = 0.035 - 0.316 (Finkelstein et al. 2009) 
and E(B - V) = 0.025 - 0.3 (Pirzkal et al. 2007). Secondly, 
since no single extinction curve (Milky Way, Small Magel- 
lanic Cloud, supernova) gives a value of f a /fc > 1; this 
larger escape fraction of Lya photons relative to the con- 
tinuum can be understood as an indication of a multi-phase 
ISM model (Neufeld 1991) where dust clumps are embedded 
in a more diffuse ionized ISM component. 

LAEs with a larger M* are also more gas rich. However, 
the ratio M g /M* is the largest (smallest) for LAEs with the 
smallest (largest) stellar mass, which implies that more mas- 
sive galaxies (Af» > 10 10 M o ) are more efficient in turning 
their gas into stars (M» 40Moyr -1 ), i.e., star forma- 
tion is suppressed (~ 8 — 25 Mq yr _1 ) in low mass galaxies 
(M» < 10 9 ' 7 M©) due to mechanical feedback. All the galax- 
ies we identify as LAEs have t« ~ 10 — 300 Myr; LAEs are 
intermediate age galaxies, instead of being extremely young 
(< 10 Myr) or extremely old (~ 1 Gyr) objects, as shown in 
Dayal et al. (2009a). The mass weighted stellar metallicity 
of LAEs scales with M» but the metallicity of the small- 
est halos has a large dispersion, possibly arising due to the 
increasing importance of feedback towards low mass halos. 
Since in our model, we have assumed SNII to be the pri- 
mary dust factories, the dust amounts are regulated solely 
by stellar processes. The dust amounts, therefore, scale with 
M», with the most star forming galaxies being most dust 
enriched. As expected, we find that f a decreases with in- 
creasing dust enrichment of the galaxy; galaxies with larger 
M*, and hence, M*, have a smaller Lya escape fraction. The 
value of f a decreases by a factor of about 3, going from 0.9 
to 0.3 as M» increases from io 85 ~ 10 - 5 M . 

We finally comment on the EW distribution obtained 
from our model and compare it to the observations (Fig. [9}. 
There is a trend of increasing EW in more luminous ob- 
jects which cannot be yet solidly identified in the avail- 
able (and uncertain) observational data. In addition, the 
observed EWs show a large dispersion, varying by as much 
as a factor of 3 for the same L a value; this suggests that the 
observed EW depends on several structural parameters of 
LAEs, namely M», peculiar velocities and dust clumping, to 
name a few. By virtue of using a LAE model which depends 
on the intrinsic galaxy properties, includes a dust calcula- 
tion and takes peculiar velocities/ density inhomogeneities 
into account, we also obtain a large spread (20-222 A) in 
the observed EWs. This is a clear improvement on our own 
previous work (Dayal et al. 2009) in which some of these 
ingredients were not yet accounted for. 

At the very last, we discuss a few caveats. First, we have 
used a constant value of the escape fraction of H i ionizing 
photons, f eac = 0.02 in all our calculations. However, the 
value of / esc remains poorly constrained both observation- 
ally and theoretically. While z ~ 3 galaxies have been used 
to obtain values of f eac < 0.04 (Fernandez-Soto et al. 2003) 
to f eac ~ 0.1 (Steidel et al. 2001), observations in the 
local Universe range between f esc < 0.01 (Deharveng et 
al. 1997) to f e3C ~ 0.1 - 0.73 (Catellanos et al. 2002). A 
larger (smaller) value of / esc would lead to larger (smaller) 



2.5 



9 2 



1.5 





42.5 43 43.5 

Log (L^/erg s~ l ) 

Figure 9. Observed EWs from Shimasaku et al. (2006) (circles) 
and model values of the observed EWs (astrexes) as a function of 
the observed Lyo luminosity. The theoretical model includes the 
full RT calculation, density, velocity fields, dust and has (xhi) ~ 
3.4 X 10~ 3 , corresponding to model M7, Tab.fl] 



H ii regions, thereby affecting the progress of the reioniza- 
tion process and hence T a ; a change in f esc would also affect 
U£\ Sec.O 

Secondly, the average age of the galaxies are calculated 
as t* = M*/M*. Since we calculate the ages assuming con- 
stant M, , it is entirely possible that some of the LAEs could 
be slightly older (younger) if M, in the past was smaller 
(larger) than the final value at z ~ 6.1. Galaxies younger 
than 10 Myr would have higher EWs as compared to the 
values shown here; however, this could be possible only for 
the the smallest galaxies, which could assemble due to merg- 
ers and undergo star formation in a time as short as 10 Myr. 

Third, in the dust model explained in Sec. 13.21 we 
have considered dust destruction by forward sweeping SNII 
shocks. However, Bianchi & Schneider (2007) have shown 
that reverse shocks from the ISM can also lead to dust de- 
struction, with only about 7% of the dust mass surviving the 
reverse shock, for an ISM density of 10~ 25 gm cm -3 . Since we 
neglect this effect, we might be over-predicting the dust en- 
richment in LAEs. This, however, does not affect our results 
because the dust optical depth depends on the surface den- 
sity of the dust distribution as shown in Eq.[4] a large (small) 
dust mass can be distributed in a large (small) volume to 
obtain identical values of the optical depth. 
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